Free thermal convection driven by nonlocal effects.

We report and explain a convective phenomenon observed in molecular dynamics simulations that cannot be classified either as a hydrodynamics instability nor as a macroscopically forced convection. Two complementary arguments show that the velocity field by a thermalizing wall is proportional to the ratio between the heat flux and the pressure. This prediction is quantitatively corroborated by our simulations. 47.27.Te, 44.25.+f Typeset using REVTEX 1 Free thermal convection —driven by buoyancy or by surface tension— is a perfectly well understood phenomenon derivable from Navier Stokes equations [1,2]. Simulations of free thermal convection by means of Molecular Dynamics (MD) techniques can be achieved with systems with as few as 10 particles and already these small systems exhibit hydrodynamic behavior as seen for example in [3][7]. Moreover MD is useful to study fluid phenomena at the microscopic level without having to make assumptions concealed behind the NavierStokes equations such as the Fourier law, Newton’s law of viscosity and local thermodynamic equilibrium. In the following we are going to present a unusual convective phenomenon (not predictable by Navier-Stokes equations) related to the variation of the temperature field in one mean free path l through the adimensional parameter l∇T/T . This parameter can be interpreted as a measure of how far from local thermal equilibrium the system is at a given point. When effects violating the local thermodynamic equilibrium are present one should question the very concept of temperature but we will manage without engaging in such delicate matters. The convective phenomenon that we are reporting takes place when there is a temperature gradient parallel to a thermal wall. The mechanism can be sketched as follows: the particles that approach a point P of the wall come from an anisotropic distribution while the particles that hit the wall at P come back to the system with a distribution which is isotropic, or at least less anisotropic than the incoming flux. A careful assessment of the difference between the incoming and outgoing fluxes at P yields the conclusion that there is a net mass flux parallel to the wall. We have observed this phenomenon in MD simulations and have made a theoretical estimation of its value. In real experiments the effect will be small but it should be observable in a rarefied gas. We have made MD simulations of a two dimensional gas of hard disks in a square box, using our own efficient algorithm [8] and the carefully devised measurement routines described in [9]. Each numerical experience consisted of two runs: (1) The system with periodic vertical walls was subjected to a temperature difference, relaxed for 200 thermal 2 diffusion times tT , and then the temperature profile T (y) was carefully measured for another 200 tT . (2) A second simulation was run under the same conditions as in (1), except that a periodic permeable thermal vertical wall was added. This new wall was defined to have the previously obtained profile T (y), namely each particle hitting the vertical wall emerged on the other side of the box with a velocity taken from a heat bath at the local temperature T (y). The sign of the vertical component of this velocity was random, and therefore microscopically it is a non slip boundary condition, in the sense that the emerging particles do not remember the velocity with which they came. Letting the particles pass through the (periodic) wall is totally irrelevant to the resulting phenomenon, but it helps reducing the boundary effects near the wall. Again the system was relaxed for 200 tT , and then measurements were averaged in time during the next 600 tT . The measurements where done dividing the system in square cells. Densities and the velocity field were measured in every cell, and fluxes were measured across the cell walls. Units are chosen such that particle’s mass and diameter, the Boltzmann constant, and the temperature at the bottom, m, D, kB and Tb respectively, are fixed to unity. With this particular choice of units the lengths are in diameter units, the temperature in energy units, and the time in units of √ mD/kBTb. The control parameters of each simulation are the number of particles N , the bulk number density nB, and the temperature at the top Tt where Tt<Tb. Our main simulation considered a system of N =1444 hard disks, bulk number density nB =0.05, implying a box side of 170 and a mean free path of about 7, and at the top the temperature was fixed to be Tt=0.1. The main observation is the following: a convective current stabilizes in the neighborhood of the vertical wall, moving towards the warmer zone. In figures 1 and 2 it is possible to see the velocity field ~v and the mass flux mn~v. At the bottom the convective current necessarily bends towards the center to come up along the central part of the box. Since the gas is highly compressible, the eye of the convective rolls are far from the expanded hotter zone. The velocity component vy in the cells by the vertical wall is almost constant, and its average 3 was vy = −0.015± 0.002 [observed] (1) after excluding 10 cells in the upper and lower extremes, with 76× 76 the total number of cells. The same convective phenomenon was observed in all the other situations we simulated: (i) n = 0.01, N = 8100, Tb = 1.0, and Tt = 0.01, (ii) n = 0.25, N = 1444, Tb = 1.0, and Tt = 0.1. The velocity component vy measured near the vertical walls in (i) was vy = 0.014±0.003 and it shows the same behavior as the preceding simulation. The vertical component of the velocity in (ii) however is not longer constant, it increases with height. The theoretical derivations that we make below are not applicable to this denser case but it is interesting to observe that the phenomenon still exists. Finally we made another simulation in which the temperature profile of the thermalizing wall was not the one obtained from a first run but rather T (y) was chosen arbitrarily to be a smooth monotonic profile. In this case we used n = 0.05, N = 1444 and Tt = 0.1. Again a convective current was created with similar characteristics to the previous ones. This result and the theoretical calculations below suggest that to obtain this convective motion, it is enough to have a temperature gradient parallel to a thermalizing wall so that each particle emerging from the wall comes from a (at least partially) thermalized distribution. This kind of convective motion can not be obtained within the usual frame of NavierStokes hydrodynamics, unless these equations are solved imposing by hand a nonvanishing tangential velocity as boundary condition at the vertical wall. That however would be artificial because hydrodynamics is built using only the first three momenta of the distribution function which are not enough information to take into consideration the phenomena that we are reporting. In what follows, we give two heuristic and complementary derivations for a rarefied two dimensional hard disk gas, one based on local nonequilibrium distribution functions, and the other one based on the mean free path theory of transport. Both derivations yield 4 essentially the same prediction for the velocity field near the vertical thermalizing wall. What the calculations below imply is that the vertical wall exerts an effective tangential force on the gas such that a velocity field — proportional to the ratio between the heat flux and the pressure — pointing against the heat flux is established. The basic idea behind the following two derivations is that particles hitting the thermalizing wall at a point P (see figure 3) come from an anisotropic nonequilibrium environment, while the particles emerging from P come from an equilibrium isotropic distribution. It is understandable then that some fluxes do not necessarily cancel and in particular we prove that there is a net velocity field. Nonequilibrium Interpretation: The velocity distribution function near a point y of the thermalizing wall (figure 3) has two contributions: (a) one from the particles that come towards y from a nonequilibrium velocity distribution and (b) the other one from the outgoing particles that come from the thermal bath at y. The nonequilibrium distribution function for a system under a heat flux adapted from [10] to the case of a two dimensional system is fneq = (